# -*- coding: utf-8 -*-
"""
Created on Sun Sep 18 09:35:03 2016

@author: zhiqiang
"""

gene_info_path = r"D:\2016试题\B\B题附件\gene_info"
name_end = '.dat'
name_begin = 'gene_'
gene_info = {}
for i in range(1,301):
    gen_info_file_Name = name_begin+str(i)+name_end
    f = open(gene_info_path+'\\'+gen_info_file_Name)
    gen_info_Data = f.readlines()
    temp = []
    for gene_wei in gen_info_Data:
        gene_wei=gene_wei.replace('\n','')
        temp.append(gene_wei)
    gene_info.setdefault(i,temp)
    f.close()
#huoqu jiyinwei xinxi

def getGeneRating(ress):
    ans = {}
    for i in range(1,301):
        ans.setdefault(i,{'bujiankang':{},'jiankang':{}})
        for weidian in gene_info[i]:
            for key in ress[weidian]['bujiankang']:
                ans[i]['bujiankang'].setdefault(key,0)
                ans[i]['bujiankang'][key] += ress[weidian]['bujiankang'][key]
            for key in ress[weidian]['jiankang']:
                ans[i]['jiankang'].setdefault(key,0)
                ans[i]['jiankang'][key] += ress[weidian]['jiankang'][key]
    return ans

gene_info_stat = getGeneRating(ress)#获取统计量
#kafang_param = {1:7.879,2:10.597,3:12.838,4:14.860,5:16.750,6:18.548,7:20.278,8:21.955,9:23.589,10:25.188,11:26.757,12:28.299,13:29.819,14:31.319,15:32.801}
kafang_param = {1: 3.841,
 2: 5.991,
 3: 7.815,
 4: 9.488,
 5: 11.071,
 6: 12.592,
 7: 14.067,
 8: 15.507,
 9: 16.919,
 10: 18.307,
 11: 19.675,
 12: 21.026,
 13: 22.362,
 14: 23.685,
 15: 24.996}
def getKaFang_3(ress,kafang=kafang_param):
    #返回卡芳分布大于的基因位
    reslist = []
    for key in ress.keys():
        res_kafang = 0
        a=[]
        b=[]
        for k in ress[key]['bujiankang'].keys():
            a.append(ress[key]['bujiankang'][k])
            if k not in ress[key]['jiankang'].keys():
                b.append(0)
            else:
                b.append(ress[key]['jiankang'][k])
        zongshu = sum(a)+sum(b)
        for i in range(len(a)):
            res_kafang += (zongshu*((a[i]-sum(a)*(a[i]+b[i])/zongshu)**2))/(sum(a)*(a[i]+b[i]))+(zongshu*((b[i]-sum(b)*(a[i]+b[i])/zongshu)**2))/(sum(b)*(a[i]+b[i]))
        if res_kafang>kafang[len(a)-1]:
#            if key=='rs12133956':
#                print(res_kafang)
            reslist.append(key)
    return reslist
third_ans = getKaFang_3(gene_info_stat)

def getWeiDian(num_list):
    temps = []
    for num in num_list:    
        temps.extend(gene_info[num])
    return temps
    
third_ans_weidian_list=getWeiDian([217,293])
#temp={}
#for key in gene_info_stat.keys():
#    temp.setdefault(str(key)+'_bujiankang',len(gene_info_stat[key]['bujiankang'].keys()))
#    temp.setdefault(str(key)+'_jiankang',len(gene_info_stat[key]['jiankang'].keys()))
#    if len(gene_info_stat[key]['bujiankang'].keys())!=len(gene_info_stat[key]['jiankang'].keys()):
#        print(key,'')
def transformPrefs(prefs):
    '''{'bujiankang': {'CC': 281, 'CT': 179, 'TT': 40},
        'jiankang': {'CC': 333, 'CT': 145, 'TT': 22}}
        转化为：{'CC': {'bujiankang': 281, 'jiankang': 333},
             'CT': {'bujiankang': 179, 'jiankang': 145},
     'TT': {'bujiankang': 40, 'jiankang': 22}}'''
    #此函数转换位点统计量
    result = {}
    for p in prefs:
        for it in prefs[p]:
            result.setdefault(it,{})
            result[it][p] = prefs[p][it]
    return result

def getBatesP(ress,weiDianList=third_ans):
    ans = {}
    for keys in weiDianList:#取出特征位点
        tranans = transformPrefs(ress[keys])
        ans.setdefault(keys,tranans)
    return ans

batesP = getBatesP(ress,third_ans_weidian_list)

def getBatesAns(person,prefs,result):
    p1=1/2
    p2=1/2
    for key in prefs.keys():#获取所有的特征key
        ks = result[key]
        k = ks[person-1]#取出person中在key这个位点上的表征
        p1 = p1*(prefs[key][k]['bujiankang']/500)
        p2 = p2*(prefs[key][k]['jiankang']/500)
    if p1>p2:
        p=1
    else:
        p=0
    return(p)


jk_bjk = [getBatesAns(i,batesP,result) for i in range(1,1001)]
jk = jk_bjk[0:500]
bjk = jk_bjk[500:1000]
print('朴素贝叶斯分类预测健康人群准确率',1-sum(jk)/500)
print('朴素贝叶斯分类预测不健康健康人群准确率',sum(bjk)/500)
print('朴素贝叶斯分类预测准确率：',(500-sum(jk)+sum(bjk))/1000)